#include "particles.h"
#include <ctime>

/*helper functions to locate particle positions and velocities
  within the X and f vectors*/
unsigned P(unsigned p, unsigned c) { return 3*p*2+c; }
unsigned V(unsigned p, unsigned c) { return 3*p*2+3+c; }

//Types of particles
enum
{
	FIREWORKS, FOUNTAIN, CHAIN, CLOTH
};

/*here starts the particle definitions*/

Particle::Particle()
{
	x = v = f = n =  vl_zero;
	m = .018;
}

Particle::Particle(Vec3d pos, Vec3d vel, float ms)
{
	x = pos;
	v = vel;
	f = vl_zero;
	m = ms;
}

void Particle::clearForces()
{
	//f = Vec3d(0,0,0);
}

/*here starts the particle system definitions*/

Vec3d addGravity(float g, float m)
{
	return Vec3d(0, g*m, 0);
}

Vec3d addSpring(Vec3d p1, Vec3d p2, float k, float r)
{
	Vec3d d = p1-p2;
	return -k*(len(d)-r)*d/len(d);
}

Vec3d addDamping(Vec3d x, float k)
{
	return -k*x;
}

ParticleSystem::ParticleSystem(unsigned num)
{
	Particle tmp;
	for (int i=0; i<num; i++) 
		particles.push_back(tmp);

	for (int i=0; i<num; i++)
		particles[i].m = 2/(float)num;
		
	f.SetSize(6*num);
	X.SetSize(6*num);
	pnum = num;
	t=0;
	ks = 30;
	vd = 1;
}

void ParticleSystem::computeAccelerations(Vecd Y, Vecd *out)
{
	static Vecd W;
	Vec3d tmp = vl_zero;
	Vec3d p1, p2;
	bool computed = false;
	float sq = sqrt(2);
	W.SetSize(6*pnum);
	for(int i=0; i<particles.size(); i++)
	{
		computed = false;
		switch(psType)
		{
		case CHAIN:
			if(i!=0)
			{
				p1 = Vec3d(Y[P(i,0)], Y[P(i,1)], Y[P(i,2)]);
				p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
				tmp = addSpring(p1, p2, ks, r);
				//p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
				//tmp += addSpring(p1, p2, ks, r);
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
				tmp += addGravity(particles[i].m, -90.81);
			}
			if(i==0)
			{
				tmp = vl_zero;
				Y[V(i,0)]=0;
				Y[V(i,1)]=0;
				Y[V(i,2)]=0;
			}
			break;
		case CLOTH:
			if(isHanged)
			if((i==0)||(i==w-1))
			{
				Y[V(i,0)]=0;
				Y[V(i,1)]=0;
				Y[V(i,2)]=0;
			}
			tmp = vl_zero;
			p1 = Vec3d(Y[P(i,0)], Y[P(i,1)], Y[P(i,2)]);
			
			//springs for the first line of particles
			if((i>=0)&&(i<=w-1))
			{
				if(i==0)
				{				
					p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
					tmp += addSpring(p2, p1, ks, r);
				}
				else
				{
					if(i==w-1)
					{
						p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						if(i+w<w*h)
						{
							p2 = Vec3d(Y[P(i+w,0)], Y[P(i+w,1)], Y[P(i+w,2)]);
							tmp += addSpring(p1, p2, ks, r);
						}
					}
					else
					{
				
						p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						if(i+w<w*h)
						{
							p2 = Vec3d(Y[P(i+w,0)], Y[P(i+w,1)], Y[P(i+w,2)]);
							tmp += addSpring(p1, p2, ks, r);
						}
					}
				}
				if(i+2<w)
				{
					p2 = Vec3d(Y[P(i+2,0)], Y[P(i+2,1)], Y[P(i+2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if(i-2>=0)
				{
					p2 = Vec3d(Y[P(i-2,0)], Y[P(i-2,1)], Y[P(i-2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
				computed = true;
			}
			//springs for the leftmost particles
			if((i%w==0)&&(i!=0)&&(!computed))
			{
				p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
				tmp += addSpring(p1, p2, ks, r);
				p2 = Vec3d(Y[P(i-w,0)], Y[P(i-w,1)], Y[P(i-w,2)]);
				tmp += addSpring(p1, p2, ks, r);
				if(i+w<w*h)
				{
					//Shear springs
					p2 = Vec3d(Y[P(i-w+1,0)], Y[P(i-w+1,1)], Y[P(i-w+1,2)]);
					tmp += addSpring(p1, p2, ks, sq*r);
				}
				if(i+2<w)
				{
					p2 = Vec3d(Y[P(i+2,0)], Y[P(i+2,1)], Y[P(i+2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if(i-2*w>=0)
				{
					p2 = Vec3d(Y[P(i-2*w,0)], Y[P(i-2*w,1)], Y[P(i-2*w,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
				computed = true;
			}
			//springs for the right most particles
			if(((i+1)%w==0)&&(i!=w-1)&&(!computed))
			{
				p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
				tmp += addSpring(p1, p2, ks, r);
				p2 = Vec3d(Y[P(i-w,0)], Y[P(i-w,1)], Y[P(i-w,2)]);
				tmp += addSpring(p1, p2, ks, r);
				if(i+w<w*h)
				{
					p2 = Vec3d(Y[P(i-w-1,0)], Y[P(i-w-1,1)], Y[P(i-w-1,2)]);
					tmp += addSpring(p1, p2, ks, sq*r);
				}
				if(w>1)
				{
					p2 = Vec3d(Y[P(i-2,0)], Y[P(i-2,1)], Y[P(i-2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if(i-2*w>0)
				{
					p2 = Vec3d(Y[P(i-2*w,0)], Y[P(i-2*w,1)], Y[P(i-2*w,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
				computed = true;
			}
			
			//springs for the last line of particles
			if((i>=w*h-w)&&(i<=w*h-1)&&(!computed))
			{
				p1 = Vec3d(Y[P(i,0)], Y[P(i,1)], Y[P(i,2)]);
				if(i==w*h-w)
				{
					p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
					tmp += addSpring(p1, p2, ks, r);
					p2 = Vec3d(Y[P(i-w,0)], Y[P(i-w,1)], Y[P(i-w,2)]);
					tmp += addSpring(p1, p2, ks, r);
					p2 = Vec3d(Y[P(i-w+1,0)], Y[P(i-w+1,1)], Y[P(i-w+1,2)]);
					tmp += addSpring(p1, p2, ks, sq*r);
				}
				else
				{
					if(i==w*h-1)
					{
						p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i-w,0)], Y[P(i-w,1)], Y[P(i-w,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i-w-1,0)], Y[P(i-w-1,1)], Y[P(i-w-1,2)]);
						tmp += addSpring(p1, p2, ks, sq*r);
					}
					else
					{
						p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i-w,0)], Y[P(i-w,1)], Y[P(i-w,2)]);
						tmp += addSpring(p1, p2, ks, r);
						p2 = Vec3d(Y[P(i-w-1,0)], Y[P(i-w-1,1)], Y[P(i-w-1,2)]);
						tmp += addSpring(p1, p2, ks, sq*r);
						p2 = Vec3d(Y[P(i-w+1,0)], Y[P(i-w+1,1)], Y[P(i-w+1,2)]);
						tmp += addSpring(p1, p2, ks, sq*r);
					}
				}
				if(i+2<w)
				{
					p2 = Vec3d(Y[P(i+2,0)], Y[P(i+2,1)], Y[P(i+2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if(i-2>=w*h-w)
				{
					p2 = Vec3d(Y[P(i-2,0)], Y[P(i-2,1)], Y[P(i-2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if(h>1)
				{
					p2 = Vec3d(Y[P(i-2*w,0)], Y[P(i-2*w,1)], Y[P(i-2*w,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
				computed = true;
				 
			}
			if(!computed)
			{
				p2 = Vec3d(Y[P(i-1,0)], Y[P(i-1,1)], Y[P(i-1,2)]);
				tmp += addSpring(p1, p2, ks, r);
				p2 = Vec3d(Y[P(i+1,0)], Y[P(i+1,1)], Y[P(i+1,2)]);
				tmp += addSpring(p1, p2, ks, r);
				p2 = Vec3d(Y[P(i-w-1,0)], Y[P(i-w-1,1)], Y[P(i-w-1,2)]);
				tmp += addSpring(p1, p2, ks, sq*r);
				p2 = Vec3d(Y[P(i-w+1,0)], Y[P(i-w+1,1)], Y[P(i-w+1,2)]);
				tmp += addSpring(p1, p2, ks, sq*r);
				if(i-2*w>0)
				{
					p2 = Vec3d(Y[P(i-2*w,0)], Y[P(i-2*w,1)], Y[P(i-2*w,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if((i-1)%w!=0)
				{
					p2 = Vec3d(Y[P(i-2,0)], Y[P(i-2,1)], Y[P(i-2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				if((i+2)%w!=0)
				{
					p2 = Vec3d(Y[P(i+2,0)], Y[P(i+2,1)], Y[P(i+2,2)]);
					tmp += addSpring(p1, p2, ks, 2*r);
				}
				tmp += addDamping(Vec3d(Y[V(i,0)], Y[V(i,1)], Y[V(i,2)]), vd);
			}
			tmp += addGravity(particles[i].m, -1000);
			tmp += particles[i].f;
			particles[i].f = vl_zero;
			break;
		case FOUNTAIN:
			tmp = addGravity(particles[i].m, -9.81);
			break;
		case FIREWORKS:
			tmp = addGravity(particles[i].m, -3.81);
			break;
		}
		W[P(i, 0)]= Y[V(i,0)];
		W[P(i, 1)]= Y[V(i,1)];
		W[P(i, 2)]= Y[V(i,2)];
		W[V(i, 0)]= tmp[0]/particles[i].m;
		W[V(i, 1)]= tmp[1]/particles[i].m;
		W[V(i, 2)]= tmp[2]/particles[i].m;
	}
	*out= W;
}

void ParticleSystem::getState()
{
	for (int i=0; i<particles.size(); i++) 
	{
		X[P(i,0)]= particles[i].x[0];
		X[P(i,1)]= particles[i].x[1];
		X[P(i,2)]= particles[i].x[2];
		X[V(i,0)]= particles[i].v[0];
		X[V(i,1)]= particles[i].v[1];
		X[V(i,2)]= particles[i].v[2];
	}
}

void ParticleSystem::setState()
{
	for (int i=0; i<particles.size(); i++) 
	{
		particles[i].t += t;
		particles[i].x[0]= X[P(i,0)];
		particles[i].x[1]= X[P(i,1)];
		particles[i].x[2]= X[P(i,2)];
		particles[i].v[0]= X[V(i,0)];
		particles[i].v[1]= X[V(i,1)];
		particles[i].v[2]= X[V(i,2)];
	}
}

void ParticleSystem::eulerStep(float dt)
{
	getState();
	computeAccelerations(X, &f);
	X = X + dt*f;
	t += dt;
	setState();
	
}

void ParticleSystem::rkStep(float dt)
{
	static Vecd k1, k2, k3, k4;
	k1.SetSize(6*pnum);
	k2.SetSize(6*pnum);
	k3.SetSize(6*pnum);
	k4.SetSize(6*pnum);
	getState();
	computeAccelerations(X, &k1);
	k1*=dt;
	computeAccelerations(X+k2/2, &k2);
	k2*=dt;
	computeAccelerations(X+k2/2, &k3);
	k3*=dt;
	computeAccelerations(X+k3, &k4);
	k4*=dt;
	X = X + k1/6 + k2/3 + k3/3 + k4/6;
	t=dt;
	setState();
	switch(psType)
	{
	case FOUNTAIN:
		for (unsigned i=0; i<pnum; i++) 
		{
			if(particles[i].x[1] < pLife)
			{
				particles[i].x = generator;
				particles[i].t = 0;
				particles[i].v[0] = 10*double(rand())/RAND_MAX - 5;
				particles[i].v[1] = double(rand())/RAND_MAX * 10.0 + 4.0;
				particles[i].v[2] = 10*double(rand())/RAND_MAX - 5;
			}
		}
		break;
	case FIREWORKS:
		for (unsigned i=0; i<pnum; i++) 
		{
			if(particles[i].t > pLife)
			{
				particles[i].x = generator;
				particles[i].t = 0;
				particles[i].v = vl_zero;
			}
		}
		break;
	}
}

void ParticleSystem::initFountain(Vec3d gen, float tt)
{
	srand(time(NULL));
	for(unsigned i=0; i<pnum; i++)
	{
		particles[i].x = gen;
		particles[i].v[0] = 10*double(rand())/RAND_MAX - 5;
		particles[i].v[1] = double(rand())/RAND_MAX * 10.0 + 4.0;
		particles[i].v[2] = 10*double(rand())/RAND_MAX - 5;
	}
	psType = FOUNTAIN;
	generator = gen;
	pLife = tt;
}

void ParticleSystem::initFireworks(Vec3d gen, float tt)
{
	srand(time(NULL));
	for(unsigned i=0; i<pnum; i++)
	{
		particles[i].x = gen;
		particles[i].t = double(rand())/RAND_MAX;
		particles[i].v[0] = 10*double(rand())/RAND_MAX-5;
		particles[i].v[1] = 10*double(rand())/RAND_MAX-5;
		particles[i].v[2] = 10*double(rand())/RAND_MAX-5;
	}
	psType = FIREWORKS;
	generator = gen;
	pLife = tt;
}

void ParticleSystem::initChain(Vec3d gen)
{
	for(int j=0; j<pnum; j++)
	{
			particles[j].x[0] = gen[0]+j;
			particles[j].x[1] = gen[1];
			particles[j].x[2] = gen[2];	
	}
	r = 1;
	generator = gen;
	psType = CHAIN;
}

void ParticleSystem::initCloth(int wt, int ht, bool hanged)
{
	int k=0;
	for(int i=0; i<ht; i++)
	{
		for(int j=0; j<wt; j++)
		{
			particles[k].x[0] = j-(float)wt/2;
			particles[k].x[1] = 5;
			particles[k].x[2] = i-(float)ht/2;
			particles[k].n = Vec3d(0,0,1);
			k++;
		}
	}
	cout << k << endl;
	isHanged = hanged;
	w = wt;
	r = 1;
	h = ht;
	psType = CLOTH;
}

void ParticleSystem::draw(bool wireframe)
{
	
	GLfloat particleColor[4];
	Vec3d tmp;
	GLfloat l;
	float wi[] = {1,1,1,1};
	switch(psType)
	{
	case FOUNTAIN:	
		particleColor[0] = 0.2;
		particleColor[1] = 0.2; 
		particleColor[2] = 1.0; 
		particleColor[3] = 1;
		glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, particleColor);
		glBegin(GL_LINES);
		
		for(int j=0; j<5; ++j)
		{
			l = j/len(particles[j].v)*.1;
				for(int i=0; i<pnum; i++)
				{
					glVertex3f(particles[i].x[0]+particles[i].v[0]*l, particles[i].x[1]+particles[i].v[1]*l, particles[i].x[2]+particles[i].v[2]*l);
					glVertex3f(particles[i].x[0]+particles[i].v[0]*.1, particles[i].x[1]+particles[i].v[1]*.1, particles[i].x[2]+particles[i].v[2]*.1);
				}
			if(j==0)
				glAccum(GL_LOAD, 1.0f/5);
			else
				glAccum(GL_ACCUM, 1.0f/5);
			glAccum(GL_RETURN, 1);
		}
		glEnd();
		break;
	case FIREWORKS:	
		particleColor[0] = 0.2;
		particleColor[1] = 0.2; 
		particleColor[2] = 1.0; 
		particleColor[3] = 1;
		
		glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, particleColor);
		glMaterialfv(GL_FRONT_AND_BACK, GL_EMISSION, wi);
		for(int j=0; j<5; ++j)
		{
			l = j/len(particles[j].v)*.1;
			glBegin(GL_LINES);
				for(int i=0; i<pnum; i++)
				{
					glVertex3f(particles[i].x[0]+particles[i].v[0]*l, particles[i].x[1]+particles[i].v[1]*l, particles[i].x[2]+particles[i].v[2]*l);
					glVertex3f(particles[i].x[0]+particles[i].v[0]*.1, particles[i].x[1]+particles[i].v[1]*.1, particles[i].x[2]+particles[i].v[2]*.1);
				}
			glEnd();
			if(j==0)
				glAccum(GL_LOAD, 1.0f/5);
			else
				glAccum(GL_ACCUM, 1.0f/5);
			glAccum(GL_RETURN, 1);
		}
		break;
	case CHAIN:
		particleColor[0] = 0.2;
		particleColor[1] = 0.2; 
		particleColor[2] = 1.0; 
		particleColor[3] = 1;
		glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, particleColor);
		for(int i=0; i<pnum; i++)
		{
			if(i!=0)
			{
				glBegin(GL_LINES);
					glVertex3f(particles[i].x[0], particles[i].x[1], particles[i].x[2]);
					glVertex3f(tmp[0], tmp[1], tmp[2]);
				glEnd();
			}
			glPushMatrix();
				glTranslatef(particles[i].x[0], particles[i].x[1], particles[i].x[2]);
				glutSolidSphere(.2, 5, 5);
			glPopMatrix();
			tmp = particles[i].x;
		}
		break;
	case CLOTH:
		particleColor[0] = 0.2;
		particleColor[1] = 0.2; 
		particleColor[2] = 1.0; 
		particleColor[3] = 1;
		float p[] = {1,1,1,1};
		glMaterialfv(GL_FRONT_AND_BACK, GL_AMBIENT_AND_DIFFUSE, particleColor);
		glMaterialfv(GL_FRONT_AND_BACK, GL_SPECULAR, p);
		glMaterialf(GL_FRONT_AND_BACK, GL_SHININESS, 128);
		int k = 0;
		glPushAttrib(GL_ALL_ATTRIB_BITS);
		if(wireframe)
		{
			glPolygonMode(GL_FRONT_AND_BACK, GL_LINE);
			glDisable(GL_LIGHTING);
			glColor3f(1,1,1);
		}
		//glBindTexture(GL_TEXTURE_2D, texture);
		glBegin(GL_TRIANGLES);
		for(int j=0; j<h-1 ; j++)
		{
			for(int i=0; i<w-1 ; i++)
			{
				glTexCoord2f(i/w, j/h);
				glNormal3f(particles[k].n[0], particles[k].n[1], particles[k].n[2]);
				glVertex3f(particles[k].x[0], particles[k].x[1], particles[k].x[2]);
				glTexCoord2f((i+1)/w, j/h);
				glNormal3f(particles[k+1].n[0], particles[k+1].n[1], particles[k+1].n[2]);
				glVertex3f(particles[k+1].x[0], particles[k+1].x[1], particles[k+1].x[2]);
				glTexCoord2f(i/w, (j+1)/h);
				glNormal3f(particles[k+w].n[0], particles[k+w].n[1], particles[k+w].n[2]);
				glVertex3f(particles[k+w].x[0], particles[k+w].x[1], particles[k+w].x[2]);
				glTexCoord2f((i+1)/w, j/h);
				glNormal3f(particles[k+1].n[0], particles[k+1].n[1], particles[k+1].n[2]);
				glVertex3f(particles[k+1].x[0], particles[k+1].x[1], particles[k+1].x[2]);
				glTexCoord2f(i/w, (j+1)/h);
				glNormal3f(particles[k+w].n[0], particles[k+w].n[1], particles[k+w].n[2]);
				glVertex3f(particles[k+w].x[0], particles[k+w].x[1], particles[k+w].x[2]);
				glTexCoord2f((i+1)/w, (j+1)/h);
				glNormal3f(particles[k+w+1].n[0], particles[k+w+1].n[1], particles[k+w+1].n[2]);
				glVertex3f(particles[k+w+1].x[0], particles[k+w+1].x[1], particles[k+w+1].x[2]);
				k++;
			} 
			k++;
		}
		glEnd();
		glPopAttrib();
		break;
	}
}

////////////////////Collition System of the Particle System //////////////////////

void ParticleSystem::computeCollition(Vec3d c, float r)
{
	for(int i=0; i<pnum; i++)
	{
		Vec3d tmp = particles[i].x-c;
		if(len(tmp)<r)
		{
			Vec3d R = c+tmp*(1+(r-len(tmp))*.1);
			Vec3d n = R/len(R);
			Vec3d tmp = Vec3d(0, dot(n, R)/n[1], 0);
			Vec3d dir = (tmp-R)/len(tmp-R);
			particles[i].x = R;
			particles[i].f = 100*dir;		
		}
			//particles[i].v *= 0;
	}
}

void ParticleSystem::computeCollition(Vec3d p, Vec3d n)
{
	for(int i=0; i<pnum; i++)
	{
		float d = dot(particles[i].x, n) - dot(n,p);
		//Vec3d tmp = particles[i].x-c;
		if(d<0)
		{
			//Vec3d R = c+tmp*(1+(r-len(tmp))*.1);
			//Vec3d n = R/len(R);
			//Vec3d tmp = Vec3d(0, dot(n, R)/n[1], 0);
			//Vec3d dir = (tmp-R)/len(tmp-R);
			particles[i].x = particles[i].x+n*(-d);
			particles[i].v = vl_zero;
			particles[i].f -= n;
			//particles[i].f = 100*dir;		
		}
			//particles[i].v *= 0;
	}
}

////////////////////////Texture part///////////////////////////////////////

struct Image {
    unsigned long sizeX;
    unsigned long sizeY;
    char *data;
};
typedef struct Image Image;


int ImageLoad(char *filename, Image *image) {
    FILE *file;
    unsigned long size;                 // size of the image in bytes.
    unsigned long i;                    // standard counter.
    unsigned short int planes;          // number of planes in image (must be 1) 
    unsigned short int bpp;             // number of bits per pixel (must be 24)
    char temp;                          // temporary color storage for bgr-rgb conversion.

    // make sure the file is there.
    if ((file = fopen(filename, "rb"))==NULL)
    {
	printf("File Not Found : %s\n",filename);
	return 0;
    }
    
    // seek through the bmp header, up to the width/height:
    fseek(file, 18, SEEK_CUR);

    // read the width
    if ((i = fread(&image->sizeX, 4, 1, file)) != 1) {
	printf("Error reading width from %s.\n", filename);
	return 0;
    }
    printf("Width of %s: %lu\n", filename, image->sizeX);
    
    // read the height 
    if ((i = fread(&image->sizeY, 4, 1, file)) != 1) {
	printf("Error reading height from %s.\n", filename);
	return 0;
    }
    printf("Height of %s: %lu\n", filename, image->sizeY);
    
    // calculate the size (assuming 24 bits or 3 bytes per pixel).
    size = image->sizeX * image->sizeY * 3;

    // read the planes
    if ((fread(&planes, 2, 1, file)) != 1) {
	printf("Error reading planes from %s.\n", filename);
	return 0;
    }
    if (planes != 1) {
	printf("Planes from %s is not 1: %u\n", filename, planes);
	return 0;
    }

    // read the bpp
    if ((i = fread(&bpp, 2, 1, file)) != 1) {
	printf("Error reading bpp from %s.\n", filename);
	return 0;
    }
    if (bpp != 24) {
	printf("Bpp from %s is not 24: %u\n", filename, bpp);
	return 0;
    }
	
    // seek past the rest of the bitmap header.
    fseek(file, 24, SEEK_CUR);

    // read the data. 
    image->data = (char *) malloc(size);
    if (image->data == NULL) {
	printf("Error allocating memory for color-corrected image data");
	return 0;	
    }

    if ((i = fread(image->data, size, 1, file)) != 1) {
	printf("Error reading image data from %s.\n", filename);
	return 0;
    }

    for (i=0;i<size;i+=3) { // reverse all of the colors. (bgr -> rgb)
	temp = image->data[i];
	image->data[i] = image->data[i+2];
	image->data[i+2] = temp;
    }
    
    // we're done.
    return 1;
}


void ParticleSystem::loadTexture(char *filename)
{
	Image *image1;
    GLuint tmp[1];
    // allocate space for texture
    image1 = (Image *) malloc(sizeof(Image));
    if (image1 == NULL) {
	printf("Error allocating space for image");
	exit(0);
    }
	
    if (!ImageLoad(filename, image1)) {
	exit(1);
    }        
	
    // Create Texture	
    /*glGenTextures(1, &tmp[0]);
    glBindTexture(GL_TEXTURE_2D, tmp[0]);   // 2d texture (x and y size)
	
    glTexParameteri(GL_TEXTURE_2D,GL_TEXTURE_MAG_FILTER,GL_LINEAR); // scale linearly when image bigger than texture
    glTexParameteri(GL_TEXTURE_2D,GL_TEXTURE_MIN_FILTER,GL_LINEAR); // scale linearly when image smalled than texture

    // 2d texture, level of detail 0 (normal), 3 components (red, green, blue), x size from image, y size from image, 
    // border 0 (normal), rgb color data, unsigned byte data, and finally the data itself.
    glTexImage2D(GL_TEXTURE_2D, 0, 3, image1->sizeX, image1->sizeY, 0, GL_RGB, GL_UNSIGNED_BYTE, image1->data);
    texture[0] = tmp[0];*/
}
